	use "${clean_data}/panel_physicians_clean.dta", clear
	
	keep if year==2017	
	gen agebin = 5*floor(age/5)
	
	preserve 
		gcollapse (mean) ptotinc income_acs totw2comp pftloss profinc wkh_imputed schedule_c_flag ///
		(count) N_UR_ptotinc = ptotinc N_UR_income_acs = income_acs N_UR_totw2comp = totw2comp N_UR_pftloss = pftloss N_UR_profinc = profinc N_UR_wkh_imputed = wkh_imputed N_UR_schedule_c_flag = schedule_c_flag , by(agebin)
		g restriction = "All ages; year=2017"
		tempfile sample1
		save `sample1'
	restore 
	
	preserve 
		keep if observed_in_acs==1
		gcollapse (mean) ptotinc income_acs (count) N_UR_ptotinc = ptotinc N_UR_income_acs = income_acs, by(agebin)
		g restriction = "All ages; year=2017; observed in ACS"
		tempfile sample2
		save `sample2'
	restore 
	
	keep if age > 54
	gcollapse (mean) received_1099ssa (count) N_UR_received_1099ssa = received_1099ssa , by(age)
	ren age agebin 
	g restriction = "age >54; year=2017"
	append using `sample1'
	append using `sample2'

	drbest_docinc ptotinc income_acs totw2comp pftloss profinc wkh_imputed schedule_c_flag received_1099ssa
	order N_UR*, last
	ds agebin restriction N_UR*, not
	foreach var of varl `r(varlist)' {
		replace `var' = . if N_UR_`var' < 4  // remove estimates with N<3 for national level per disclosure rules
	}
	export delimited using "${mypath}/intermediate_csv/desc03-lifecycle_earnings_binscatter.csv", replace dataf delim(tab)
